Introduction

We wish to create conservation-sensitive, climate-smart fisheries closures for the ABNJ of the Pacific Ocean using prioritizr.

Prioritizr (Hanson et al., 2018) is a spatial prioritization software (R package) that uses integer linear programming (ILP) allowing it to create exact solutions to large spatial planning problems faster than its heuristic counterparts (e.g. Marxan).

This file runs multiple scripts found in the main Github repository: tinbuenafe/PacificProj, where the codes are explained in greater detail.

Assembling the data layers

Before running prioritizr, we first had to assemble the following data layers:

  1. Study Area
  2. Fisheries conservation features
  3. Climate smart features
  4. Cost layer

Study Area

The study area has the following characteristics:

  • Pacific ABNJ
  • Limited boundaries to areal jurisdiction of RFMOs responsible for managing large pelagic fisheries (IATTC and WCPFC)
  • Robinson projection

We ran the following code using polygons from rnaturalearth package and the EEZ .shp files from Marine Regions (Flanders Marine Institute, 2019) to create a Pacific-centered map of the global ABNJs.

The relevant scripts for creating the study area are labeled with a prefix of 01_StudyArea in ~/scripts/.

# to create the Pacific-centered maps:
source("scripts/01a_StudyArea_CreatingPUs.R")
# this script uses the following functions (embedded actual script):

# function to create polygon of the ABNJ by masking the EEZs and the landmasses
source("scripts/study_area/fCreateMaskedPolygon.R") 
# function to reproject these polygon to Robinson projection
source("scripts/study_area/fConvert2PacificRobinson.R") 
# function to create the hexagonal planning units
source("scripts/study_area/fCreate_PlanningUnits.R") 

This figure shows the Pacific-centered global map (using Robinson projection).

world_robinson <- readRDS("outputs/01_StudyArea/01a_StudyArea/PacificCenterLand.rds")
global_plot <- ggplot() +
  geom_sf(data = world_robinson, colour = "grey20", fill="grey20", size = 0.1, show.legend = "line") +
  theme_bw()
global_plot

We then created equal-area hexagonal planning units, or PUs, (0.5\(^\circ\) x 0.5\(^\circ\) at the equator) of the study area. This resulted in 31,917 PUs.

# define objects for plotting; objects are the result of running 01_StudyArea
pacific_robinson <- readRDS("outputs/01_StudyArea/01a_StudyArea/PacificCenterABNJ.rds")
PUsPac <- readRDS("outputs/01_StudyArea/01a_StudyArea/PacificABNJGrid_05deg.rds")
Bndry <- fCreateRobinsonBoundary(west = 78, east = 140, north = 51, south = 60) # coordinates in degrees

study_area <- ggplot() +
    geom_sf(data = pacific_robinson, colour = NA, fill = NA, size = 0.2, show.legend = "line") +
    geom_sf(data = world_robinson, color = "grey20", fill="grey20", size = 0.1, show.legend = "line") +
    geom_sf(data = PUsPac, colour = "grey64", aes(fill = "ABNJ"), size = 0.1, show.legend = TRUE) + 
    scale_fill_manual(name = "Study Area", values = c("ABNJ" = "coral3")) +
    coord_sf(xlim = c(st_bbox(Bndry)$xmin, st_bbox(Bndry)$xmax), # Set limits based on Bndry bbox
            ylim = c(st_bbox(Bndry)$ymin, st_bbox(Bndry)$ymax),
            expand = TRUE) +
    theme_bw()
study_area

Fisheries conservation features

The following were used in this study:

  1. Distribution of selected bycatch species
  2. Distribution of larval spawning grounds of selected commercially targeted species

Bycatch distribution

First, we compiled a list of bycatch reported in the Pacific ABNJ. Then, we extracted their IUCN distribution. These polygons are from the IUCN Red List, showing the known ranges of the species (IUCN Red List, 2021).

We limited the bycatch included in this study to sea turtles, but this workflow can be easily replicated and expanded to include more species.

The following were the sea turtles reported as bycatch and intersected with the study area:

  • Loggerhead turtle (Caretta caretta)
  • Green turtle (Chelonia mydas)
  • Leatherback turtle (Dermochelys coriacea)
  • Hawksbill turtle (Eretmochelys imbricata)
  • Olive ridley turtle (Lepidochelys olivacea)

The relevant scripts for creating the bycatch distribution layer are labeled with a prefix of 04_IUCN in ~/scripts/.

The plots below show the distribution of these species using fIUCNINtersect().

# the function is found in:
source("scripts/04a_IUCN_fIUCNINtersect.R")
# the function is ran in:
source("scripts/04b_IUCN_fIUCNIntersectRun.R")

Spawning grounds of commercially targeted species

Just like the bycatch layer, we compiled a list of commercially targeted large pelagic species managed in the Pacific from different sources. From this list, we have 4 species with available data on their global larval spawning grounds:

  • Yellowfin Tuna (Thunnus albacares)
  • Albacore (Thunnus alalunga)
  • Skipjack Tuna (Katsuwonus pelamis)
  • Swordfish (Xiphias gladius)

The models for these species are from a report by Mercer (2019) where he created Generalized Additive Models (GAMs) to predict spawning grounds. He used average seasonal abundance from a historical dataset of longline fisheries surveys (dated back to the 1950s), and oceanographic and environmental conditions from global databases.

We reran his models (see scripts labeled with a prefix 05a_Commercial in ~/scripts/ for complete reruns) and used the medians of the predictions for each species as probability thresholds to create the distribution maps and the function fCommercialFeat().

# to intersect the results of the GAMs with the PUs
# the function is in:
source("scripts/05b_Commercial_fCommercialFeat.R")
# the function is ran at:
source("scripts/05c_Commercial_fCommercialFeatRun.R")
global_plots

This is an example of how the data looks like (so far) for bycatch:

head(bycatch)
## Simple feature collection with 6 features and 2 fields
## Geometry type: POLYGON
## Dimension:     XY
## Bounding box:  xmin: -3779074 ymin: 168708.9 xmax: -3695824 ymax: 473118.3
## CRS:           +proj=robin +lon_0=180 +x_0=0 +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m +no_defs
## # A tibble: 6 x 3
##   cellsID                                                geometry species       
##     <int>                                           <POLYGON [m]> <chr>         
## 1       1 ((-3751324 216773.5, -3779074 232795.1, -3779074 26483… Caretta_caret…
## 2       2 ((-3751324 312902.8, -3779074 328924.4, -3779074 36096… Caretta_caret…
## 3       3 ((-3751324 409032.1, -3779074 425053.6, -3779074 45709… Caretta_caret…
## 4       4 ((-3723574 168708.9, -3751324 184730.5, -3751324 21677… Caretta_caret…
## 5       5 ((-3723574 264838.2, -3751324 280859.7, -3751324 31290… Caretta_caret…
## 6       6 ((-3723574 360967.5, -3751324 376989, -3751324 409032.… Caretta_caret…

And for the large pelagics:

head(pelagic)
## Simple feature collection with 6 features and 2 fields
## Geometry type: POLYGON
## Dimension:     XY
## Bounding box:  xmin: -3779074 ymin: 168708.9 xmax: -3695824 ymax: 473118.3
## CRS:           +proj=robin +lon_0=180 +x_0=0 +y_0=0 +ellps=WGS84 +datum=WGS84 +units=m +no_defs
## # A tibble: 6 x 3
##   cellsID                                                       geometry species
##     <int>                                                  <POLYGON [m]> <chr>  
## 1       1 ((-3751324 216773.5, -3779074 232795.1, -3779074 264838.2, -3… ALB    
## 2       2 ((-3751324 312902.8, -3779074 328924.4, -3779074 360967.5, -3… ALB    
## 3       3 ((-3751324 409032.1, -3779074 425053.6, -3779074 457096.7, -3… ALB    
## 4       4 ((-3723574 168708.9, -3751324 184730.5, -3751324 216773.5, -3… ALB    
## 5       5 ((-3723574 264838.2, -3751324 280859.7, -3751324 312902.8, -3… ALB    
## 6       6 ((-3723574 360967.5, -3751324 376989, -3751324 409032.1, -372… ALB

Each sf object has the following columns:

  • cellsID is the unique ID for each hexagonal PU
  • species is the species code (e.g. ALB = Albacore)
  • geometry

Climate-smart features

To create climate-smart closures, we prioritized protection of areas with:

  1. low levels of exposure to climate warming, and
  2. high levels of potential retention of biodiversity.

We determined these areas using the following metrics:

  1. Relative Climate Exposure (RCE) index

    \(\large RCE(yr^{-1}) = \frac{{\text Slope \hspace{2mm}} (^{\circ}Cyr^{-1} \hspace{2mm} 2050 - 2100)}{{\text Seasonal \hspace{2mm} Range} \hspace{2mm} (^{\circ}C \hspace{2mm} 2015 - 2020)}\)

  • numerator is the slope of the linear regression of projected mean annual temperatures
  • denominator is the current mean seasonal temperature range

Areas with low RCE index values are areas that will be exposed to low levels of warming relative to their current seasonal fluctuations of temperature. Areas of low RCE index values are prioritized for selection.

  1. Climate velocity

    \(\large Climate \hspace{2mm} velocity \hspace{2mm} (kmyr^{-1}) = \frac{{\text Slope \hspace{2mm}} (^{\circ}Cyr^{-1} \hspace{2mm} 2050 - 2100)}{{\text Spatial \hspace{2mm} gradient} \hspace{2mm} (^{\circ}Ckm^{-1} \hspace{2mm} 2050 - 2100)}\)

  • numerator is the slope of the linear regression of projected mean annual temperatures
  • denominator is calculated from the vector sum of the latitudinal and longitudinal pairwise differences of the mean temperature at each focal cell using a 3 x 3 neighborhood window

Areas of slow climate velocity are areas that are more likely to retain their current environmental conditions and consequentially, their biodiversity. Areas of slow climate velocity are prioritized for selection.

The temperatures used to calculate the aforementioned metrics were derived from 12 Global Circulation Models (GCMs) from the Couple Model Intercomparison Project Phase 6 (CMIP6; https://esgf-node.llnl.gov). The GCMs are under these three climate scenarios:

  1. SSP1-2.6: optimistic scenario
  2. SSP2-4.5: intermediate scenario
  3. SSP5-8.5: pessimistic scenario

RCE

We intersected the RCE values with the PUs using fClimateInt(). The relevant scripts for creating the climate-smart features layer are found in ~/scripts/ with the prefix 07_Climate.

# function is found at:
source("scripts/07a_Climate_fClimateInt.R")
# function is ran at:
source("scripts/07b_Climate_fClimateIntRun.R")

Climate velocity

We intersected the climate velocity values with the PUs using the same function, fClimateInt().

Cost layer

The cost layer represents the opportunity cost to fishing of establishing the proposed fisheries closures.

The cost layer for each planning unit was represented in US Dollars, and was calculated by multiplying the:

  • estimated total catch of medium - large pelagics in kg (from Watson, 2017),and
  • the ex-vessel price of the same species in the aforementioned functional groups in US$ \(kg^{-1}\) (from Tai et al., 2017)
# function is found at:
source("scripts/06a_Cost_fCostPU.R")
# function is ran at:
source("scripts/06b_Cost_fCostPURun.R")
## [1] "minimum: 0.0036497397813946 maximum: 5977.41455078125"

We used log-transformed cost for the figure below:

Representation targets

Filtering and retaining PUs

We join the climate-smart layer with the conservation features layer. Then, we only retain planning units that have RCE and/or climate-velocity values belonging to their lower 25th percentile range. This was how we priortized protection of areas that are predicted to have lower exposure to climate warming and higher retention of biodiversity.

We used the function fFilterQuartile() and the scripts that are relevant for this step are found in ~/scripts/ with the prefix 08_Filter.

# function is found at:
source("scripts/08a_Filter_fFilterQuartile.R")
# function is ran at:
source("scripts/08d_Filter_fFilterQuartile25Run_NoProv.R")

This was done for both the bycatch and large pelagic layers, across different scenarios (SSP1-2.6, SSP2-4.5, SSP5-8.5). Below is an example of how the data looks like, with the following columns:

  • cellsID: unique ID for each PU
  • species: species code
  • total_area: the total area covered by each subtyped feature (in km\(^{-1}\))
  • velocity: climate velocity
  • velo_tvalue: transformed values of climate velocity (velocity * 10)
  • RCE: RCE index values
  • RCE_tvalue: transformed values of RCE (cube root of RCE)
# for bycatch under the climate scenario SSP1-2.6.
bycatchSSP126 <- readRDS("outputs/08_Filter/08d_Filter25_NoProv/bycatchIUCNSSP126_25percentile.rds") %>% 
  as_tibble() %>% 
  dplyr::select(-area_km2, -cellsID.2, -geometry)
head(bycatchSSP126)
## # A tibble: 6 x 7
##   cellsID new_features        velocity velo_tvalue     RCE RCE_tvalue total_area
##     <int> <chr>                  <dbl>       <dbl>   <dbl>      <dbl>      <dbl>
## 1      12 Caretta_caretta_IU…   0.239        2.39  0.0352       0.328  19761581.
## 2      21 Caretta_caretta_IU…   0.143        1.43  0.0239       0.288  19761581.
## 3      22 Caretta_caretta_IU…   0.0163       0.163 0.00493      0.170  19761581.
## 4      23 Caretta_caretta_IU…  -0.0532      -0.532 0.0107       0.221  19761581.
## 5      32 Caretta_caretta_IU…   0.237        2.37  0.0355       0.329  19761581.
## 6      33 Caretta_caretta_IU…   0.0717       0.717 0.0131       0.236  19761581.
# for large pelagics under the cliamte scenario SSP5-8.5.
commercialSSP585 <- readRDS("outputs/08_Filter/08d_Filter25_NoProv/commercialSSP585_25percentile.rds") %>% 
  as_tibble() %>% 
  dplyr::select(-area_km2, -cellsID.2, -geometry)
head(commercialSSP585)
## # A tibble: 6 x 7
##   cellsID new_features velocity velo_tvalue   RCE RCE_tvalue total_area
##     <int> <chr>           <dbl>       <dbl> <dbl>      <dbl>      <dbl>
## 1      81 ALB              35.8        358.  2.29       1.32   4278830.
## 2      84 ALB              37.1        371.  2.24       1.31   4278830.
## 3      85 ALB              52.9        529.  2.34       1.33   4278830.
## 4      86 ALB              38.7        387.  1.18       1.06   4278830.
## 5      88 ALB              39.1        391.  2.19       1.30   4278830.
## 6      89 ALB              46.5        465.  1.29       1.09   4278830.

Aside from these climate-smart closures, we also explored climate-uninformed, wherein the climate-smart features were not included in the analyses. Subsequently, the step of filtering and retention of planning units was also not included.

source("scripts/08e_Filter_fFilterQuartile100Run_NoProv.R")

This is how the data looks like:

# for bycatch:
bycatch_uninformed <- readRDS("outputs/08_Filter/08e_Filter100_NoProv/bycatchIUCNuninformed_100percentile.rds") %>% 
  as_tibble() %>% 
  dplyr::select(-area_km2, -geometry)
head(bycatch_uninformed)
## # A tibble: 6 x 3
##   cellsID new_features         total_area
##     <int> <chr>                     <dbl>
## 1       1 Caretta_caretta_IUCN  62747287.
## 2       2 Caretta_caretta_IUCN  62747287.
## 3       3 Caretta_caretta_IUCN  62747287.
## 4       4 Caretta_caretta_IUCN  62747287.
## 5       5 Caretta_caretta_IUCN  62747287.
## 6       6 Caretta_caretta_IUCN  62747287.
# for large pelagics:
commercial_uninformed <- readRDS("outputs/08_Filter/08e_Filter100_NoProv/commercialuninformed_100percentile.rds") %>%
  as_tibble() %>% 
  dplyr::select(-area_km2, -geometry)
head(commercial_uninformed)
## # A tibble: 6 x 3
##   cellsID new_features total_area
##     <int> <chr>             <dbl>
## 1       1 ALB            9582019.
## 2       2 ALB            9582019.
## 3       3 ALB            9582019.
## 4       4 ALB            9582019.
## 5       5 ALB            9582019.
## 6       6 ALB            9582019.

Assigning targets

To determine the relative targets per subtyped feature, we checked each one’s IUCN category (IUCN, 2015) using rredlist(). Threatened features (i.e. classified as vulnerable, endangered, threatened) were assigned a fixed maximum target. Features (\(i\)) that were not classified as threatened (e.g. classified as near-threatened, etc.) were assigned a relative target based on their distribution.

Relative target\(_i\) \(= Target_{max} \cdot \frac{PU_i}{PU_{total}} \cdot (Target_{max} - Target_{min})\)

  • PU_i: planning units of feature \(i\)
  • PU_{total}: total planning units
  • Target_{max}: Maximum target / the fixed maximum target assigned to threatened features (e.g. 1)
  • Target_{min}: Minimum target (e.g. 0.1)

Lower targets are assigned to features with broader distribution, and higher targets to features with more restricted distribution.

We reran this step for an array of targets, for both climate-smart and climate-uninformed closures:

Runs Maximum Minimum
1 1 0.1
2 0.9 0.1
3 0.8 0.1
4 0.7 0.1
5 0.6 0.1
6 0.5 0.1
7 0.4 0.1
8 0.3 0.1
9 0.2 0.1
10 0.1 0

We used the function fRepresentTarget(), and the scripts are found in ~/scripts/ with the prefix 09.

# the function is found in:
source("scripts/09a_Target_fRepresentTarget.R")
# the runs for climate-smart closures are found in:
source("scripts/09d_Target_fRepresentTargetSmartRun_NoProv.R")
# the runs for climate-uninformed closures are found in:
source("scripts/09e_Target_fRepresentTargetUninformedRun_NoProv.R")

For this .rmd, we’ll be showing results using Runs 2 (Max = 0.9, Min = 0.1), 5 (Max = 0.6, Min = 0.1) and 8 (Max = 0.3, Min = 0.1). This workflow shows that it can be replicated using other parameters.

Below is an example of how the data looks like:

target90_commercial <- readRDS("outputs/09_Target/09d-e_TargetRuns/02_Target90/SSP126/target_commercialSSP126.rds") %>% 
  as_tibble()
head(target90_commercial)
## # A tibble: 4 x 5
##   new_features total_area category                               geometry target
##   <chr>             <dbl> <chr>                        <MULTIPOLYGON [m]>  <dbl>
## 1 ALB            3603928. NT       (((-2169566 -3772591, -2197317 -37565…   3.05
## 2 SKP            2238116. LC       (((-2086316 -2859363, -2114066 -28433…   1.89
## 3 SWO            4390870. LC       (((-1448063 -2907428, -1475813 -28914…   3.71
## 4 YFT            3942713. NT       (((-948560.7 -696454.6, -976310.8 -68…   3.33
target90_bycatch <- readRDS("outputs/09_Target/09d-e_TargetRuns/02_Target90/SSP126/target_bycatchSSP126.rds") %>% 
  as_tibble()
head(target90_bycatch)
## # A tibble: 5 x 5
##   new_features total_area category                               geometry target
##   <chr>             <dbl> <chr>                        <MULTIPOLYGON [m]>  <dbl>
## 1 Rep-1347         56020. VU       (((-1753314 -2859363, -1781065 -28433…     90
## 2 Rep-2226         53352  EN       (((-1753314 -2859363, -1781065 -28433…     90
## 3 Rep-3437      20644556. VU       (((688697.1 -5743241, 660946.9 -57272…     90
## 4 Rep-4127         58687. CR       (((-1753314 -2859363, -1781065 -28433…     90
## 5 Rep-5414         16006. VU       (((-1864315 -2763234, -1892065 -27472…     90
  • new_features: features
  • total_area: the total area covered by each subtyped feature (in km\(^{-1}\))
  • category: IUCN redlist category
  • target: assigned target in percentage

Running prioritizr

prioritizr (v. 7.0.1.2) uses the following functions to define the problem.

# installing latest version of prioritizr from GitHub
devtools::install_github("prioritizr/prioritizr")
# installing gurobi as the optimizer
install.packages("/Library/gurobi911/mac64/R/gurobi_9.1-1_R_4.0.2.tgz", repos = NULL)
library(priortizr)
library(gurobi)

p1 <- prioritizr::problem(x, feature_list, "cost") %>% 
    add_min_set_objective() %>%
    add_relative_targets(target_list) %>%
    add_binary_decisions() %>%
    add_gurobi_solver(gap = 0.1, verbose = FALSE) 

s1 <- solve(p1)
  • problem(): we define the data
  • x: data frame having all the data layers (all planning units with their specific ID (cellsID), features data, cost)
  • feature_list: list of the features in x
  • cost: name of the column in x where the cost is found
  • add_min_set_objective(): function that allows us to employ the minimum set objective function. This allows prioritizr to create spatial plans that will minimize the cost.
  • add_relative_targets(): target_list list of assigned targets for each feature
  • add_binary_decisions(): prioritizr will only either select or not select a PU
  • add_gurobi_solver(): Gurobi was used as the optimizer to solve the problem (at a default 10%)
  • solve(): prioritizr solves the problem p1

Prioritizr runs

We used the fPrioritizrRun() function to streamline the runs. fPrioritizrRun() arranges all the data and makes sure that they’re in the right format, then runs problem() and solve(). This function also saves the summaries each plan (cost, % area protected and number of representative features that achieved 30% or more protection).

The relevant scripts for this step is found in ~/scripts/ with the prefix 10_Prioritizr.

# the function is in:
source("scrpits/10a_Prioritizr_fPrioritizr.R")
# to run the function with climate-smart data
source("scripts/10h_Prioritizr_fPrioritizrIUCNSmart_NoProv.R")
# to run the function with climate-uninformed data
source("scripts/10i_Prioritizr_fPrioritizrIUCNUninformed_NoProv.R")
# to create no-regret closures with the climate-smart plans

To create no-regret plans, we used thre function fCreateNoRegret() which intersected the selected PUs in all three climate scenarios. The relevant scripts for this step is found in ~/scripts with the prefix 11_NoRegret.

# the function is in:
source("scripts/11_NoRegret_fCreateNoRegret.R")
# the runs are found in:
source("scripts/11e_NoRegret_fCreateNoRegretIUCNRuns_NoProv.R")

As mentioned in the previous section, we will only be showing results of Runs 2 (Max = 0.9, Min = 0.1), 5 (Max = 0.6, Min = 0.1) and 8 (Max = 0.3, Min = 0.1).

Run 2: Maximum target = 90%

Climate-smart runs

SSP1-2.6 run

SSP126_target90 <- readRDS("outputs/10_Prioritizr/10h-i_IUCNRuns_NoProv/02_Target90/solution_SSP126.rds")

  • cost is in USD
  • perc_area_protected is the percentage of the Pacific ABNJ protected
  • protected_PU is the number of features that have >= 30% (by area) protection
## # A tibble: 0 x 12
## # … with 12 variables: target_scenarios <chr>, sum_area <dbl>,
## #   percent_area <dbl>, num_pu <dbl>, total_cost <dbl>, median_velocity <dbl>,
## #   median_tvelocity <dbl>, median_RCE <dbl>, median_RCEtvalue <dbl>,
## #   represented_features <dbl>, target <dbl>, scenario <chr>

SSP2-4.5 run

SSP245_target90 <- readRDS("outputs/10_Prioritizr/10h-i_IUCNRuns_NoProv/02_Target90/solution_SSP245.rds")

  • cost is in USD
  • perc_area_protected is the percentage of the Pacific ABNJ protected
  • protected_PU is the number of features that have >= 30% (by area) protection
## # A tibble: 0 x 12
## # … with 12 variables: target_scenarios <chr>, sum_area <dbl>,
## #   percent_area <dbl>, num_pu <dbl>, total_cost <dbl>, median_velocity <dbl>,
## #   median_tvelocity <dbl>, median_RCE <dbl>, median_RCEtvalue <dbl>,
## #   represented_features <dbl>, target <dbl>, scenario <chr>

SSP5-8.5 run

SSP585_target90 <- readRDS("outputs/10_Prioritizr/10h-i_IUCNRuns_NoProv/02_Target90/solution_SSP585.rds")

  • cost is in USD
  • perc_area_protected is the percentage of the Pacific ABNJ protected
  • protected_PU is the number of features that have >= 30% (by area) protection
## # A tibble: 0 x 12
## # … with 12 variables: target_scenarios <chr>, sum_area <dbl>,
## #   percent_area <dbl>, num_pu <dbl>, total_cost <dbl>, median_velocity <dbl>,
## #   median_tvelocity <dbl>, median_RCE <dbl>, median_RCEtvalue <dbl>,
## #   represented_features <dbl>, target <dbl>, scenario <chr>

No-regret plan

noregret_target90 <- readRDS("outputs/11_NoRegret/11e_IUCNRuns_NoProv/noregretclosures_Target90.rds")

Climate-uninformed run

uninformed_target90 <- readRDS("outputs/10_Prioritizr/10h-i_IUCNRuns_NoProv/02_Target90/solution_uninformed.rds")

  • cost is in USD
  • perc_area_protected is the percentage of the Pacific ABNJ protected
  • protected_PU is the number of features that have >= 30% (by area) protection
## # A tibble: 1 x 12
##   target_scenarios       sum_area percent_area num_pu total_cost median_velocity
##   <chr>                     <dbl>        <dbl>  <dbl>      <dbl>           <dbl>
## 1 Target90_summaryuninf…   6.72e7         79.0  25207    767309.       67242193.
## # … with 6 more variables: median_tvelocity <dbl>, median_RCE <dbl>,
## #   median_RCEtvalue <dbl>, represented_features <dbl>, target <dbl>,
## #   scenario <chr>

Run 5: Maximum target = 60%

Climate-smart runs

SSP1-2.6 run

SSP126_target60 <- readRDS("outputs/10_Prioritizr/10h-i_IUCNRuns_NoProv/05_Target60/solution_SSP126.rds")

  • cost is in USD
  • perc_area_protected is the percentage of the Pacific ABNJ protected
  • protected_PU is the number of features that have >= 30% (by area) protection
## # A tibble: 0 x 12
## # … with 12 variables: target_scenarios <chr>, sum_area <dbl>,
## #   percent_area <dbl>, num_pu <dbl>, total_cost <dbl>, median_velocity <dbl>,
## #   median_tvelocity <dbl>, median_RCE <dbl>, median_RCEtvalue <dbl>,
## #   represented_features <dbl>, target <dbl>, scenario <chr>

SSP2-4.5 run

SSP245_target60 <- readRDS("outputs/10_Prioritizr/10h-i_IUCNRuns_NoProv/05_Target60/solution_SSP245.rds")

  • cost is in USD
  • perc_area_protected is the percentage of the Pacific ABNJ protected
  • protected_PU is the number of features that have >= 30% (by area) protection
## # A tibble: 0 x 12
## # … with 12 variables: target_scenarios <chr>, sum_area <dbl>,
## #   percent_area <dbl>, num_pu <dbl>, total_cost <dbl>, median_velocity <dbl>,
## #   median_tvelocity <dbl>, median_RCE <dbl>, median_RCEtvalue <dbl>,
## #   represented_features <dbl>, target <dbl>, scenario <chr>

SSP5-8.5 run

SSP585_target60 <- readRDS("outputs/10_Prioritizr/10h-i_IUCNRuns_NoProv/05_Target60/solution_SSP585.rds")

  • cost is in USD
  • perc_area_protected is the percentage of the Pacific ABNJ protected
  • protected_PU is the number of features that have >= 30% (by area) protection
## # A tibble: 0 x 12
## # … with 12 variables: target_scenarios <chr>, sum_area <dbl>,
## #   percent_area <dbl>, num_pu <dbl>, total_cost <dbl>, median_velocity <dbl>,
## #   median_tvelocity <dbl>, median_RCE <dbl>, median_RCEtvalue <dbl>,
## #   represented_features <dbl>, target <dbl>, scenario <chr>

No-regret plan

noregret_target60 <- readRDS("outputs/11_NoRegret/11e_IUCNRuns_NoProv/noregretclosures_Target60.rds")

Climate-uninformed run

uninformed_target60 <- readRDS("outputs/10_Prioritizr/10h-i_IUCNRuns_NoProv/05_Target60/solution_uninformed.rds")

  • cost is in USD
  • perc_area_protected is the percentage of the Pacific ABNJ protected
  • protected_PU is the number of features that have >= 30% (by area) protection
## # A tibble: 1 x 12
##   target_scenarios       sum_area percent_area num_pu total_cost median_velocity
##   <chr>                     <dbl>        <dbl>  <dbl>      <dbl>           <dbl>
## 1 Target60_summaryuninf…   4.48e7         52.7  16807    393372.       44834353.
## # … with 6 more variables: median_tvelocity <dbl>, median_RCE <dbl>,
## #   median_RCEtvalue <dbl>, represented_features <dbl>, target <dbl>,
## #   scenario <chr>

Run 8: Maximum target = 30%

Climate-smart runs

SSP1-2.6 run

SSP126_target30 <- readRDS("outputs/10_Prioritizr/10h-i_IUCNRuns_NoProv/08_Target30/solution_SSP126.rds")

  • cost is in USD
  • perc_area_protected is the percentage of the Pacific ABNJ protected
  • protected_PU is the number of features that have >= 30% (by area) protection
## # A tibble: 0 x 12
## # … with 12 variables: target_scenarios <chr>, sum_area <dbl>,
## #   percent_area <dbl>, num_pu <dbl>, total_cost <dbl>, median_velocity <dbl>,
## #   median_tvelocity <dbl>, median_RCE <dbl>, median_RCEtvalue <dbl>,
## #   represented_features <dbl>, target <dbl>, scenario <chr>

SSP2-4.5 run

SSP245_target30 <- readRDS("outputs/10_Prioritizr/10h-i_IUCNRuns_NoProv/08_Target30/solution_SSP245.rds")

  • cost is in USD
  • perc_area_protected is the percentage of the Pacific ABNJ protected
  • protected_PU is the number of features that have >= 30% (by area) protection
## # A tibble: 0 x 12
## # … with 12 variables: target_scenarios <chr>, sum_area <dbl>,
## #   percent_area <dbl>, num_pu <dbl>, total_cost <dbl>, median_velocity <dbl>,
## #   median_tvelocity <dbl>, median_RCE <dbl>, median_RCEtvalue <dbl>,
## #   represented_features <dbl>, target <dbl>, scenario <chr>

SSP5-8.5 run

SSP585_target30 <- readRDS("outputs/10_Prioritizr/10h-i_IUCNRuns_NoProv/08_Target30/solution_SSP585.rds")

  • cost is in USD
  • perc_area_protected is the percentage of the Pacific ABNJ protected
  • protected_PU is the number of features that have >= 30% (by area) protection
## # A tibble: 0 x 12
## # … with 12 variables: target_scenarios <chr>, sum_area <dbl>,
## #   percent_area <dbl>, num_pu <dbl>, total_cost <dbl>, median_velocity <dbl>,
## #   median_tvelocity <dbl>, median_RCE <dbl>, median_RCEtvalue <dbl>,
## #   represented_features <dbl>, target <dbl>, scenario <chr>

No-regret plan

noregret_target60 <- readRDS("outputs/11_NoRegret/11e_IUCNRuns_NoProv/noregretclosures_Target60.rds")

Climate-uninformed run

uninformed_target30 <- readRDS("outputs/10_Prioritizr/10h-i_IUCNRuns_NoProv/08_Target30/solution_uninformed.rds")

  • cost is in USD
  • perc_area_protected is the percentage of the Pacific ABNJ protected
  • protected_PU is the number of features that have >= 30% (by area) protection
## # A tibble: 1 x 12
##   target_scenarios       sum_area percent_area num_pu total_cost median_velocity
##   <chr>                     <dbl>        <dbl>  <dbl>      <dbl>           <dbl>
## 1 Target30_summaryuninf…   2.24e7         26.3   8407    149748.       22426513.
## # … with 6 more variables: median_tvelocity <dbl>, median_RCE <dbl>,
## #   median_RCEtvalue <dbl>, represented_features <dbl>, target <dbl>,
## #   scenario <chr>

Summary

Since we ran prioritizr from 10 - 100% maximum target representation (in increments of 10%), we’re showing a summary of how cost and percent protection varied with increasing target representation.

source("scripts/12d_Stat_fCreateSummaryRuns.R")
# showing the log(cost), log(% area protected), number of protected planning units, and number of adequately represented features.
# this was done to determine which maximum representation targets to report/use in the runs
smart_plans

We also determined Cohen’s Kappa Coefficient between plans created for each run. Below are the summaries for each run as well as a matrix plot of the correlations between plans. We also show the summaries of the representation of each feature for each spatial plan.

The relevant scripts for this part are found in ~/scripts/ with prefix 12_Stat.

# the function for creating the correlation plots is found in:
source("scripts/12e_Stat_fKappaCorrplot.R")
# the runs are found in:
source("scripts/12f_Stat_fKappaCorrplotRun.R")

# the function for creating circular bar plots to show summaries of the representation of each feature:
source("scripts/12g_Stat_fCreateCircBarPlot.R")
# the runs are found in:
source("scripts/12h_Stat_fCreateCircBarPlotRuns.R")

Summary for Run 2: Maximum Target = 90%

## # A tibble: 4 x 12
##   target_scenarios       sum_area percent_area num_pu total_cost median_velocity
##   <chr>                     <dbl>        <dbl>  <dbl>      <dbl>           <dbl>
## 1 Target90_summarySSP126   2.25e7         26.4   8428    266081.          0.0541
## 2 Target90_summarySSP245   3.00e7         35.2  11246    324104.          2.31  
## 3 Target90_summarySSP585   3.46e7         40.7  12975    354942.         13.1   
## 4 Target90_summaryuninf…   6.72e7         79.0  25207    767309.   67242193.    
## # … with 6 more variables: median_tvelocity <dbl>, median_RCE <dbl>,
## #   median_RCEtvalue <dbl>, represented_features <dbl>, target <dbl>,
## #   scenario <chr>

Summary for Run 5: Maximum Target = 60%

## # A tibble: 4 x 12
##   target_scenarios       sum_area percent_area num_pu total_cost median_velocity
##   <chr>                     <dbl>        <dbl>  <dbl>      <dbl>           <dbl>
## 1 Target60_summarySSP126   1.46e7         17.2   5480    124930.          0.0455
## 2 Target60_summarySSP245   1.83e7         21.4   6846    158994.          2.29  
## 3 Target60_summarySSP585   2.06e7         24.2   7734    165253.         12.9   
## 4 Target60_summaryuninf…   4.48e7         52.7  16807    393372.   44834353.    
## # … with 6 more variables: median_tvelocity <dbl>, median_RCE <dbl>,
## #   median_RCEtvalue <dbl>, represented_features <dbl>, target <dbl>,
## #   scenario <chr>

Summary for Run 8: Maximum Target = 30%

## # A tibble: 4 x 12
##   target_scenarios       sum_area percent_area num_pu total_cost median_velocity
##   <chr>                     <dbl>        <dbl>  <dbl>      <dbl>           <dbl>
## 1 Target30_summarySSP126   7.16e6         8.41   2684     43083.          0.0576
## 2 Target30_summarySSP245   9.50e6        11.2    3561     60609.          2.34  
## 3 Target30_summarySSP585   1.07e7        12.6    4026     62372.         13.1   
## 4 Target30_summaryuninf…   2.24e7        26.3    8407    149748.   22426513.    
## # … with 6 more variables: median_tvelocity <dbl>, median_RCE <dbl>,
## #   median_RCEtvalue <dbl>, represented_features <dbl>, target <dbl>,
## #   scenario <chr>